/*
 *  VectorHelper.h
 *  NSE_SOLVER
 *
 *  Created by Patrick Zulian on 1/15/11.
 *  Copyright 2011 ICS - USI - Lugano. All rights reserved.
 *
 */

#ifndef NSE_SOLVER_VECTOR_HELPER_H
#define NSE_SOLVER_VECTOR_HELPER_H

#include "Base.h"
#include "Utils.h"

namespace nsesolver {
	
	
	/**
	 * VectorType Interface has to implement operator[] and getSize()
	 * VectorType and T Interface has to support all arithmetic operators and its op= version
	 * T also T t = 0
	 */
	template< typename T >
	class VectorHelper {
	public:
		template< class VectorType >
		static inline 
		void random(VectorType &v, const T &range)
		{
			for (int i = 0; i < v.size(); i++) {
				v[i] = static_cast< T > ( range * utils::Rand() ); 
			}
			
		}
		
		template< class VectorType >
		static inline 
		void random(VectorType &v, const T min, const T max)
		{
			const T range = max - min;
			for (int i = 0; i < v.size(); i++) {
				v[i] =  min + static_cast< T >( range *  utils::Rand() ); 
			}
		}
		
		template< class VectorType >
		static inline 
		T dotProduct(const VectorType &v1, const VectorType &v2)
		{
			NSE_SOLVER_ASSERT( v1.size() == v2.size() );
			
			T result = 0;
			for (int i = 0; i < v1.size(); i++) {
				result += static_cast< T >( v1[i] * v2[i] );
			}
			return result;
		}
		
		template< class VectorType >
		static inline 
		T squaredNorm(const VectorType &v)
		{
			T result = 0;
			for (int i = 0; i < v.size(); i++) {
				result += static_cast< T >( v[i] * v[i] );
			}
			return result;
		}
		
		template< class VectorType >
		static inline 
		T norm(const VectorType &v)
		{
			return static_cast< T > (sqrt( squaredNorm(v) ));
		}
		
		template< class VectorType >
		static inline 
		T squaredDistance(const VectorType &v1, const VectorType &v2)
		{
			NSE_SOLVER_ASSERT( v1.size() == v2.size() );
			
			return squaredNorm( v1 - v2 );
		}
		
		
		template< class VectorType >
		static inline 
		T distance(const VectorType &v1, const VectorType &v2)
		{
			NSE_SOLVER_ASSERT( v1.size() == v2.size() );
			
			return norm( v1 - v2 );
		}
		
		
		template< class VectorType >
		static inline
		void allSet(VectorType &v, const T &value)
		{
			for (int i = 0; i < v.size(); i++) {
				v[i] = value;
			}
		}
		
		
		template< class VectorType >
		static inline
		void normalize(VectorType &v)
		{
			v /= norm(v);
		}
		
		template< class VectorType >
		static inline
		const VectorType normalized(const VectorType &v)
		{
			VectorType result(v);
			normalize(result);
			return result;
		}
		
		
		template< class VectorType >
		static inline 
		void inPlaceSum(VectorType &result, const VectorType &v2)
		{
			NSE_SOLVER_ASSERT( result.size() == v2.size() );
			
			for (int i = 0; i < result.size(); i++) {
				result[i] += v2[i];
			}

		}
		
		template< class VectorType >
		static inline 
		void inPlaceSubtract(VectorType &result, const VectorType &v2)
		{
			NSE_SOLVER_ASSERT( result.size() == v2.size() );
			
			for (int i = 0; i < result.size(); i++) {
				result[i] -= v2[i];
			}

		}
		
		
		template< class VectorType >
		static inline 
		void inPlaceMultiply(VectorType &result, const VectorType &v2)
		{
			NSE_SOLVER_ASSERT( result.size() == v2.size() );
			

			for (int i = 0; i < result.size(); i++) {
				result[i] *= v2[i];
			}
		}
		
		template< class VectorType >
		static inline 
		void inPlaceScale(VectorType &result, const T &scaleFactor)
		{
			for (int i = 0; i < result.size(); i++) {
				result[i] *= scaleFactor;
			}

		}
		
		
		//New Vector version of inPlace operations
		template< class VectorType >
		static inline 
		const VectorType sum(const VectorType &v1, const VectorType &v2)
		{
			VectorType result(v1);
			inPlaceSum(result, v2);
			return result;
		}
		
		template< class VectorType >
		static inline 
		const VectorType subtract(const VectorType &v1, const VectorType &v2)
		{
			VectorType result(v1);
			inPlaceSubtract(result, v2);
			return result;
		}
		
		template< class VectorType >
		static inline 
		const VectorType multiply(const VectorType &v1, const VectorType &v2)
		{
			VectorType result(v1);
			inPlaceMultiply(result, v2);
			return result;
		}
		
		
		template< class VectorType >
		static inline 
		const VectorType scale(const VectorType &v, const T &scaleFactor)
		{
			VectorType result(v);
			inPlaceScale(result, scaleFactor);
			return result;
		}
		
	};
	
	
}

#endif



